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ABSTRACT 



We estimate the galaxy stellar mass function and stellar mass density for star-forming and quiescent galaxies with 0.2 < z < 4. We 
construct a large, deep (K^ < 24) sample of 220,000 galaxies selected using the new UltraVISTA DRl data release. Our analysis is 
based on precise 30-band photometric redshifts. By comparing these photometric redshifts with 10,800 spectroscopic redshifts from 
the zCOSMOS bright and faint surveys, we find a precision of tr^;/(i+j) = 0.008 at i* < 22.5 and o-a-hi+-<i = 0.03 at 1.5 < z < 4. 
We derive the stellar mass function and correct for the Eddington bias. We find a mass-dependent evolution of the global and star- 
forming populations, with the low-mass end of the mass functions evolving more rapidly than the high-mass end. This mass-dependent 
evolution is a direct consequence of the star formation being "quenched" in galaxies more massive than M > 10'" AIq, while 
the rapid evolution at the low mass end is explained by the evolution of the specific star formation rate. By deriving the global stellar 
mass density, we show that galaxies grow in mass twice as quickly at 1 < z < 4 than at 0. 1 < z < 1 . We also confirm that global 
stellar mass density is lower by 0.2-0.3 dex than the cosmic star fomation rate integrated over cosmic time. For the mass function 
of the quiescent galaxies, we do not find any significant evolution of the high-mass end at z < 1, while however we observe a clear 
flattening of the faint-end slope. From z ~ 3 to z ~ 1, the density of quiescent galaxies increases over the entire mass range. Their 
comoving stellar mass density increases by 1.5 dex between z ~ 3 and z ~ 1 and by less than 0.2 dex at z < 1. Finally, we compare our 
results with the semi-analytical model of Wang et al. (2008) and find that these models overestimate the density of low mass quiescent 
galaxies by an order of magnitude, while the density of low-mass star-forming galaxies is successfully reproduced. 

Key words. Galaxies: distances and redshifts - Galaxies: evolution - Galaxies: formation - Galaxies: star formation - Galaxies: 
stellar content 



* Based on data products from observations made with ESO 

Telescopes at the La Silla Paranal Observatory under ESO programme Cambridge Astronomy Survey Unit on behalf of the UltraVISTA con- 
ID 179.A-2005 and on data products produced by TERAPIX and the sortium. 
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1. Introduction 

The galaxy stellar mass function (hereafter MF) is a fundamental 
indicator of the physical processes that regulate mass assembly 
in galaxies across cosmic time. Stellar mass assembly in galaxies 
is believed to result from several physical processes such as star 
formation from in-situ or accreted gas as well as from major 
or minor mergers. The star formation may also be stopped, or 
"quenched", by external process such as "feedback" in Active 
Galactic Nuclei (AGN) for galaxies residing in more massive 
dark matter haloes, or by other processes such as supemovae- 
driven winds in less massive haloes. The relative contribution 
and operating timescales of these different processes are still a 
matter of debate. 

The evolutionary tracks of the MF as a function of look-back 
time reveal the major paths taken by different galaxy populations 
across cosmic time. Despite the considerable progress made over 
the last decade, our understanding of the early evolution of stel- 
lar mass growth is still incomplete and several puzzling issues 
have yet to be resolved: ( i) The stellar mass density of the most 
massive galaxies (lO'^'^A!©) at z ~ 3 is almost identical to the 
local measurement (Perez-Gonzalez et al. 2008). This result re- 
quires a powerful feedback mechanism in order to halt star for- 
mation in massive galaxies, like AGN feedback (e.g. Bower et 
al. 2006, Croton et al. 2006). However, numerous systematic un- 
certainties could still be hiding an evolution of the high-mass 
end (Marchesini et al. 2009). In particular, if both random and 
systematic errors are considered, Marchesini et al. (2009) "can 
not exclude a strong evolution (by as much as a factor of ~50) in 
the number density of the most massive galaxies (> lO^'^Al©) 
from z = 4.0 to z = 1.3". (ii) The stellar mass density of the most 
massive quiescent galaxies undergoes little evolution (less than 
0.2 dex) from z ~ 1 to z ~ 0. 1 (Arnouts et al. 2007, Pozzetti et al. 
2010, Ilbert et al. 2010). Therefore, major merger activity must 
be relatively limited for this population since z ~ 1, in agreement 
with galaxy merger rate measurements (e.g. Lopez-Sanjuan et 
al., 2012). (Hi) The stellar mass density increases continuously 
with cosmic time. The star formation rate integrated along cos- 
mic time and the stellar mass density evolution should provide 
a coherent picture (Arnouts et al. 2007, Boissier et al. 2010), 
but are still difficult to reconcile unless the initial mass function 
(IMF) changes with time, as advocated by several authors (e.g. 
Wilkins et al. 2008, Lu et al. 2012). A significant uncertainty 
remains on the contribution of low mass galaxies which could 
impact the global stellar mass density at z > 2 (Santini et al. 
2012). (iv) Red galaxies are created very efficiently at 1 < z < 2 
(Cirasuolo et al. 2007, Arnouts et al. 2007, Ilbert et al. 2010, 
Cassata et al. 2011). This implies that the quenching of star 
forming galaxies must be extremely efficient at z > 1 . Arnouts et 
al. (2007), Ilbert et al. (2010) and Kajisawa et al. (2011) found 
that the stellar mass density of quiescent galaxies increased by 
roughly one order of magnitude from z ~ 2 to z ~ 1 (2.5Gyr) 
as star-formation stops in galaxies which then migrate into the 
"red sequence" of passively evolving galaxies. But the amount 
of evolution is still debated. Brammer et al. (2011) find also an 
evolution of the quiescent population at 1 < z < 2, but at a 
smaller rate of 0.5 dex for the massive population. Ilbert et al. 

(2010) identify z ~ 1 as a transition epoch in the assembly of the 
most massive quiescent galaxies (the density of the most massive 
galaxies ceases to increase at z < 1 ). But from the Brammer et al. 

(201 1) study, no clear transition occurs at z ~ 1 . Therefore, a key 
point is to consolidate these measurements at z < 2 and follow 
the growth of the quiescent population to even earlier times. 



Beyond z ~ 1 observations are challenging, and are currently 
limited either by depth or small numbers of objects in relatively 
small fields, making the computation of the massive end of the 
MF sensitive to cosmic variance and the lower mass end diffi- 
cult to constrain. The identification of a robust quiescent galaxy 
sample requires accurate photometric redshifts with a low num- 
ber of catastrophic failures. A field of at least one square degree 
or more is necessary to provide a large volume and minimise 
cosmic variance. Moreover, samples need to be deep enough to 
probe the low-mass end of the MF beyond z ~ 2. The availabil- 
ity of deep wide field near-IR multi-band photometry is essen- 
tial. Spectral features like the D4000 or Balmer break move into 
the near-IR at z > 1.5, and several near-IR bands are required 
to properly sample the spectral energy distribution (SED) and 
enable a stable photometric redshift and type classification from 
SED-fitting techniques. 

The COSMOS field is one of the best available fields to 
derive the MF thanks to the large area (2 deg^) and the large 
amount of deep (Iab ~ 26.5) multi-wavelength data available 
(more than 35 bands). Several papers showing the MF evolution 
in the COSMOS field have been published (Drory et al. 2009, 
Ilbert et al. 2010, Pozzetti et al. 2010, Dominguez-Sanchez et al. 
2011). New photometric and spectroscopic datasets have been 
obtained in the last two years which allow us to greatly reduce 
the systematic uncertainties in the MF estimate at z > 1 . The first 
Ultra VISTA DRl data released (McCracken et al. 2012) covers 
1.5 deg^ in four near-infrared filters Y, J, H and K^. The DRl 
data are at least one magnitude deeper than previous COSMOS 
near-infrared datasets (McCracken et al. 2010) and also provide 
new Y-band photometric information. Since the Balmer break 
lies between UltraVISTA filters at z > 1.3, we can now de- 
rive photometric redshifts which are more robust at high red- 
shift. Secondly, almost 35,000 new spectra are now available in 
the COSMOS field. This sample includes 9,900 spectroscopic 
redshifts at z > 1, including extremely faint objects. Such a 
spectroscopic sample is essential to ensure that our analysis can 
be extended at high redshift. New NIR spectroscopic samples 
(WFMOS and MOIRCS on Subaru and WFC3/HST grism data) 
contain spectroscopic redshifts for the quiescent and the dusty 
populations at z > 1 .5. Therefore, we can ensure that our photo-z 
are robust at z > 1 and we are now able to extend MF measure- 
ments to z = 4 and to lower stellar masses than previous studies 
in the same field. 

In this paper, we extend MF and stellar mass density mea- 
surements out to z = 4 and for stellar masses down to 10'"-^A<o 
for the global population. The new data sets are introduced in 
^ The photometric redshifts and associated physical proper- 
ties are discussed in |3] The method used to estimate the MF 
and the associated uncertainties is given in ^ We present the 
measured MF and stellar mass density for the full, star-forming 
and quiescent galaxy samples in ^ Results are discussed and 
compared to semi-analytical models in ^ In particular, we in- 
vestigate systematic uncertainties linked to the choice of the star 
formation histories in our models. Throughout this paper, we 
use the standard cosmology (0,„ - 0.3, Qa = 0.7 with 
H{) - 70 km s"' Mpc"'). Magnitudes are given in the AB sys- 
tem (Oke 1974). The stellar masses are given in units of solar 
masses (Mq) for a Chabrier (2003) IMF 



' www.eso.org/sci/observing/phase3/data_releases/ultravista_drl.htnil 
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2. Data description 

2.1. Preparation of stacked images and confidence maps 

Our photometric catalogue comprises near-infrared data taken 
with the VIRCAM (Emerson & Sutheriand 2010) on the VISTA 
telescope as part of the UltraVISTA project and optical broad 
and intermediate-band data taken with the SUPRIME camera 
on Subaru in support of the COSMOS project (Capak et al. 
2007). The near-infrared data we use here corresponds to the 
UltraVISTA DRl data release fully described in McCracken et 
al. (2012). 

To construct our multi-band catalogue, we first down- 
loaded all COSMOS tiles and confidence maps from the IRSAQ 
COSMOS archive, choosing where possible the "best" seeing 
images. Since the IRSA tiles have the same tangent point and 
pixel scale as UltraVISTA DRl, they may be simply "pasted to- 
gether" (i.e., without image resampling) using the swarp soft- 
ware to produce a single, large image which is pixel-matched 
to UltraVISTA DRl. Catalogues can then be simply generated 
using SExtractor (Bertin & Amouts 1996) in "dual-image" 
mode, using matched apertures on each image. The list of broad 
and intermediate band / narrow band images used is as follows: 
u% Bj, Vj, r+, /+, z+, /A484, /A527, M624, M679, M738, 
/A767, M427, /A464, /A505, /A574, /A709, /A827, A^B711, 
NB816. Table 1. in Ilbert et al (2010) lists the effective wave- 
lengths of each of these filters (note that we do not use the 
Subaru g-band data as this has particularly poor seeing). All 
four UltraVISTA bands are used: Y, J, H, K^; the depths for 
the UltraVISTA DRl ai-e given in Table 1 of McCracken et al. 
(2012). 

Before catalogue extraction, we must first take into account 
the large variation in seeing between different COSMOS images 
and also within the UltraVISTA stacks themselves. As already 
documented in McCracken et al. (2012), in the UltraVISTA H 
and stacks there are "columns" of different seeing as a con- 
sequence of the different VISTA pawprints contained in the final 
stacks being taken under different seeing conditions. To correct 
for this, we first construct for each H and stack six sepa- 
rate stacks comprising the six individual VISTA "pawprints". 
We then measure the average seeing on each of six pawprints us- 
ing the PSFex software (Bertin et al. 2011), which corresponds 
to a fit of a Moffat (1969) profile. Next, the individual pawprints 
are convolved by a Gaussian profile calculated to bring the PSF 
on the final images to 1". A similar procedure is adopted to ho- 
mogenise the COSMOS broad and intermediate-band images: 
the seeing is measured on each stack and the images are then 
convolved with a Gaussian profile to bring the final image to 1 " 
which is the worst PSF among all bands (the Y band from Ultra- 
VISTA). 

2.2. IRAC observations 

The IRAC data consist of all cryogenic data (Sanders et al. 2007) 
and the data from the Warm mission SEDS program (Fazio et al. 
2012) covering ~ 0. Ideg^ in the center of the field. The data were 
reduced with the MOPEXQ software package. Photometric mea- 
surements were made using SExtractor in dual image mode 
with the Subaru ;-band image as a detection image and the IRAC 
image for measurement. This improves photometric accuracy 
when optical sources are close to each other by making use of 
the SExtractor aperture de-blending routines. As in described 

^ |http : //irsa . i pac . caltech . edu/Mi ssions /cosmos . html | 
^ |http : // ssc ■ spitzer . caltech . edu/postbcd71 



in Sanders et al. (2007) photometric measurements were made 
in a 3.8"diameter aperture and corrected to pseudo-total using 
a statistical aperture correction. A comparison with the fluxes 
measured in Ilbert et al. (2010) shows only a small systematic 
offset, expected as a consequence of the update to the IRAC cal- 
ibration. 

2.3. Source catalogue extraction and merged catalogue 
creation 

In extracting our source catalogues using SEXtractor, the 
choice of the detection image is important. Our scientific ob- 
jectives drive us to use the longest possible wavelengths in order 
to reliably detect galaxies at intermediate to high redshifts. On 
the other hand, we also want a sample which is as complete as 
possible to the faintest possible limits in stellar mass. We con- 
struct a detection image using SWarp (Bertin et al. 2002) from 
the chi-squared sum of the (non-convolved) UltraVISTA DRl 
YJHKs images, following the techniques outlined in Szalay et 
al. (1999). This ensures that all sources detected in at least one 
VISTA band are included in the final catalog. 

With this detection image and a set of PSF-homogenised im- 
ages in hand, we can proceed to catalogue extraction. In order to 
make reliable magnitude error estimates on stacked, convolved 
data which has a non-negligible amount of correlated noise we 
use for each Subaru band pre-computed effective gain values 
(Capak et al. 2007) in combination with the MAP_RMS confi- 
dence maps supplied by IRSA. This is particularly important for 
the shorter wavelength data which saturates at relatively bright 
magnitudes. For the UltraVISTA data, where we use weight 
maps as opposed to RMS-maps we set the SEXtractor op- 
tion RESCALE_WEIGHTS to N to ensure that image convolution 
has no effect on image noise measurement and that the weight- 
maps are not rescaled automatically as is normally the case in 
SExtractor. To account for additional noise sources not ac- 
counted for in the RMS maps such as errors in the background 
subtraction, zero-point errors, confusion, morphological aper- 
ture corrections, and other effects the errors were multiplied by 
a factor of 1 .5 to match the measured noise on the extracted pho- 
tometry. 

For each detected source we measure aperture magnitudes 
into a 3" diameter circle, and "pseudo-total" "Kron" magnitudes 
(Kron, 1980), corresponding to SExtractor's MG_APER and 
MAG_AUTO respectively. Next, catalogues from each band were 
merged together into a single FITS table and galactic extinction 
values computed at each object position using the Schlegel et 
al. (1998) dust maps were added. In addition object mask flags 
indicating bad regions in optical and near-infrared bands were 
included, and saturated pixels in the optical bands were flagged 
by using the appropriate FLAGJlAPs at the extraction stage. 

In each band, for the purposes of photometric redshift mea- 
surements, aperture magnitudes are corrected to pseudo-total 
magnitudes using a band-dependent aperture correction com- 
puted using a stellar curve-of-growth method (although the im- 
ages are PSF-homogenised, there are still small residual PSF 
variations band-to-band). The SED fitting is performed on these 
magnitudes. However, this aperture correction is not enough to 
capture the total flux of large and bright galaxies. As a con- 
sequence, we would underestimate the stellar masses of these 
galaxies. To get the correct stellar masses, we compute an aver- 
age transformation between these aperture corrected magnitudes 
and SExtractor's MAG_AUTO over all broad-band filters for each 
object and then apply this adjustment to the stellar masses. 



3 



O. Ilbert et al.: Mass assembly in quiescent and star-forming galaxies since z - 4 from Ultra VISTA 



^ 3 



o 

I— I 



o 







1 1 1 1 

" 12482 


1 1 1 1 1 1 1 
galaxies 


1 1 




; rj: 2.0% 






- "az/(i+z) 

- 


0.0097 






_ □ 

□ 


" ■'"□°°'^J 


-■■ 






^ □ -- °^^PP^^B 




— 












b oB 


a □ 


/COSMOS fainl 






° X 


zCOSMOS bright 




^"*4lll^^ .-X □ ° 




DEIMOS 






A 


Comparat t'0RS2 - 


^ III 






^ilvcrma,n l''N40S 






o 








Onodera MOIRCS - 






A 


Krogager WFC3 














□ ^ °f ° 






1 1 1 1 1 1 1 ° 


1° ° i„ 


: ° 1° , , , , 



12 3 4 

spectroscopic redshift 

Fig. 1. Photometiic redshifts versus spectroscopic redshifts. 
Only secure spectroscopic redshifts at /Tj < 24 are considered. 
Different symbols correspond to the spectroscopic samples la- 
beled in the right. 

Finally, the /- selected IRAC catalogue described above was 
cross-matched to the near-infrared Ultra VISTA catalogue. This 
/- selected catalogue also contains GALEX UV magnitudes 
computed by Zamoski et al (2007). GALEX fluxes were mea- 
sured using PSF fitting method. 

3. Photometric redshifts and physical parameters 

3.1. Photometric redshift estimation 

We derive our photometric redshifts following the recipes out- 
lined in Ilbert et al. (2009). The photometric redshifts are de- 
rived using "Le Phare" (Amouts et al. 2002, Ilbert et al. 2006). 
Ilbert et al. (2009) used 31 templates including elliptical and 
spiral galaxies from the Polletta et al. (2006) library and 12 
templates of young and blue star-forming galaxies generated 
with Bruzual and Chariot stellar population synthesis models 
(2003, hereafter BC03). Using a spectroscopic sample of qui- 
escent galaxies, Onodera et al. (2012) showed that the estimate 
of the photo-z for the quiescent galaxies in Ilbert et al. (2009) 
were underestimated at 1.5 < z < 2. We improve the photo-z 
for this specific population by adding two new templates of el- 
liptical galaxies generated with BC03. These two new templates 
are generated assuming an exponentially declining star forma- 
tion rate with a short timescale r = Q.3Gyr and two different 
metalHcities (Z = 0.008 and Z = 0.02 i.e. Z0). They include 22 
ages well sampled between O.SGyr and 4Gyr old. Thanks to a 
better sampling of younger ages than our previous library, they 
improve the photo-z for the quiescent population at z > 1 .5. 

Extinction is added as a free parameter (E(B-V)<0.5) and 
several extinction laws are considered (Calzetti et al. 2000, 
Prevot 1984 and a modified version of the Calzetti laws includ- 
ing a bump at 2175A). We do not add any extinction for the 
templates redder than SO. We also do not allow additional ex- 
tinction for the templates redder than Sc (the Sa, Sb, Sc tem- 



plates from Polletta et al. 2007 already include some extinction). 
Emission lines are added to the templates using an empirical re- 
lation between the UV light and the emission line fluxes (Ilbert 
et al. 2009). By contrast with Ilbert et al. (2009), we assign the 
redshift using the median of the marginalised probability distri- 
bution function rather than the minimum While the results 
are broadly similar between the two methods, using the median 
has the advantage of producing more reliable error bars for the 
photometric redshift^ and reducing the effect of aliasing in the 
photometric redshift space. 

We combine several spectroscopic samples to test the ac- 
curacy of the photometric redshifts, including zCOSMOS- 
bright with 20700 bright VIMOSA'^LT spectra selected at /+ < 
22.5 (Lilly et al. 2007), zCOSMOS faint with 9500 faint 
VIMOS/VLT spectra selected at 1.5 < z < 3 (Lifly et al., 
in preparation), 2300 DEIMOS/Keck redshifts which combined 
several selected sub-populations of blue star-forming and in- 
frared galaxies at 0.5 < z < 6 (Kartaltepe et al. 2012, Capak et 
al., in preparation), 835 FORS2/VLT redshifts at 0.6 < z < 1.8 
(Comparat et al., in preparation), 138 FMOS/Subaru redshifts 
at 1.4 < z < 1.8 (Silverman et al., in preparation), 18 faint 
quiescent galaxies at z < 1.9 obtained with MOIRCS/Subaru 
(Onodera et al. 2012) and 16 faint quiescent galaxies at 1.85 < 
z < 2.6 obtained with the WFC3 grism observations from the 
3D-HST survey (Krogager et al., in preparation). We keep only 
the most secure spectroscopic redshifts (e.g. flag 3 and 4 for 
zCOSMOS) at < 24 which reduces to 12482 the number of 
spectroscopic redshifts used for the comparison. The compari- 
son between photometric and spectroscopic redshifts is shown in 
Figureliland the accuracy obtained for each spectroscopic sam- 
ple is listed in Table 1 . 

FigurelU shows photometric redshift accuracy has two 
regimes. At z < 1.5, the Balmer break falls between the inter- 
mediate band filters. In this regime, the spectroscopic sample is 
dominated by the zCOSMOS bright sample selected at < 22.5 
(blue crosses). The precision is better than 1 % with less than 1 % 
of catastrophic failures. Even as faint as ~ 24, comparisons 
with the DEIMOS sample show a better than 3% precision at 
z < 1.5 (green crosses). At intermediate redshifts 1.4 < z < 2, 
the FMOS sample from Silverman et al. shows an excellent pre- 
cision of 0.026. However we note that //„ is identified in the NIR 
spectra from a knowledge of photo-z which could explain the 
low failure rate. At higher redshifts 1.5 < z < 3.5, the accuracy 
of the photo-z has been tested against the zCOSMOS faint sam- 
ple. Here the precision is around 3%, showing that we can safely 
extend our analysis at z > 1.5 for faint star-forming galaxies. 
Thanks to new NIR spectroscopic samples, we can now test the 
accuracy of the photo-z for the quiescent population at z > 1.5. 
We combine the quiescent sample of Onodera et al. (2012) and 
Krogager et al. (in preparation). These galaxies are extremely 
faint with a median band magnitude of 24.5 and a median 
redshift of 1 .8. We have only one outlier in 30 galaxies. The pre- 
cision is as good as o-az/(\+z,) - 0.056. We still detect a small 
bias with medianizp - zs) - -0.1, which does not have signifi- 
cant impact on our work. 

We also compute the best-fitting using the stellar library. 
The sources are classified as stars when the obtained with the 
star library is lower than the one obtained with the galaxy library. 
Since near-infrared data are available, such criterion is effective 
in isolating the stars. 



We integrate the PDF over 68% of its area around the median solu- 
tion 
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Table 1. Characteristics of the spectroscopic redshift samples 
and photometric redshift accuracy. Only the most secure spec- 
troscopic redshifts at K^^ < 24 are considered. The median red- 
shift and magnitude are provided for each sample. The frac- 
tion of catastrophic failures t] is defined as galaxies with \Zphoi - 
Zspec\/{i + Zspec) > 0.15. The redshift precision o-^-/([+-^) is com- 
puted using the normalized median absolute deviation (Ilbert et 
al. 2006). 



3.2. Stellar masses and associated completeness 

We turn now to the estimation of physical properties based 
on these photometric redshifts, first estimating galaxy stellar 
masses. We rely on a model to convert the galaxy luminosity 
into stellar mass. We generate a library of synthetic spectra nor- 
maUsed at one solar mass. These synthetic spectra are fitted to 
the multi-colour photometry described above using LeJ'hare. 
The physical parameter called "stellar mass" in this paper cor- 
responds to the median of the marginalization over the stel- 
lar mass probability distribution. The library of synthetic spec- 
tra is generated using the Stellar Population Synthesis (SPS) 
model of Bruzual and Chariot (2003). Several other models are 
available in the literature (e.g. Fioc & Rocca-Volmerange 1997, 
Maraston 2005, Bruzual 2007, Conroy et al. 2009). The stellar 
masses could vary by 0.1-0.15 dex depending on the considered 
SPS model (e.g. Walcher et al. 2011). For consistency, we use 
the same default library as Ilbert et al. (2010): we assume the 
Calzetti (2000) extinction law; emission line contributions are 
included using an empirical relation between the UV light and 
the emission line fluxes (Ilbert et al. 2009); we used three differ- 
ent metallicities (Z = 0.004, Z = 0.008, Z = 0.02 i.e. Zo); the 
star formation history declines exponentially following r^'e"'^^. 
We considered nine possible t values ranging from 0. 1 Gyr to 
30 Gyr. We will test the systematic uncertainties linked to our 
choice of star formation history in ^6.31 Following Fontana et al. 
(2006), Pozzetti et al. (2007), Ilbert et al. (2010), we impose the 
prior E(B - V) < 0.15 if age It > 4 (a low extinction is imposed 
for galaxies having a low SFR). 

Figure |2] shows the stellar masses as a function of redshift. 
We can detect galaxies with masses as low as M ~ 10"'Mo at 
z = 4. The mass limit depends on the mass-to-light ratio of the 
considered template. Figure |2] shows the stellar mass limit that 
can be reached for different templates. For the clarity, we show 
only two extreme templates in terms of star formation history 
(t - O.lGyr and t - 3QGyr), with solar metallicity and we se- 
lect three ages. Oldest galaxies have the higher stellar mass limit. 
We follow a procedure similar to Pozzetti et al. (2010) to define 
the stellar mass completeness limit. We base our estimate on the 
90% of the templates which are the most often fitted. In practice, 
we compute the lowest stellar mass which could be detected for 
a galaxy with log{Mumit) - log{M)+Q.A-{K-2A), given a sample 
selected ?X K < 2A (brown points). At a given redshift, the stel- 
lar mass completeness limit corresponds to the mass with 90% 
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Fig. 2. Stellar masses versus redshift for the Kab < 24 selected 
sample. The top and bottom panels correspond to the full sam- 
ple and the quiescent population, respectively. Black points cor- 
respond to the measured stellar masses. The brown points cor- 
respond to Miimii which is the lowest mass which could be ob- 
served for a given galaxy for a limit at Kab - 24. The blue cir- 
cles coiTespond to the completeness limit chosen in this study. 
The red and green curves correspond to the two extreme tem- 
plates with T = O.lGyr and t - iQGyr, respectively. We show 
3 ages which are 0.9, 2 and 5 Gyr (from the bottom to the top, 
respectively). The dashed lines correspond to an extinction of 
E{B - V) - 0.2. If the lines are not drawn, it means that either 
the age is older than the age of the Universe or that the condition 
age It > 4 implies a small extinction. We also indicate in the top 
of each panel the fraction of galaxies detected in at least 4 bands. 
The fraction is given in % per redshift bin. The top and bottom 
numbers coiTespond to a detection threshold at > 2cr and > 5cr, 
respectively. 



of the galaxies having their Ainmii below the stellar mass com- 
pleteness limit. Following this procedure, not more than 10% of 
the galaxies could be missed in the lowest mass end of the mass 
function. 

At z > 2, galaxies could be too faint to be detected in opti- 
cal and we need to rely on optical upper-limits and NIR fluxes 
to estimate the photo-z. Photometric redshifts with a precision 
of 5% are routinely computed with four bands (e.g. Ilbert et al. 
2006). Therefore, we checked that most of the galaxies are de- 
tected in four bands at least, to allow a robust photo-z estimate. 
The fraction of galaxies detected in a minimum of four bands at 
2 and 5 cr are given in Figure |2] This fraction is always above 
95% at z < 2 even if we consider a 5cr detection limit. At z > 2, 
the fraction stays higher than 90% (80%) for a 2cr (5cr) detection 
limit, except at 2.5 < z < 3 for the quiescent (still 86% have their 
magnitudes measured in a minimum of 4 bands, with a detection 
limit better than 2cr). 
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3.3. Galaxy classification 

In order to divide the sample into quiescent and star-forming 
galaxies, we use slightly modified version of the two-colour se- 
lection technique proposed by Williams et al. (2009). Following 
Ilbert et al. (2010), we use the rest-frame two-colour selection 
NUV - versus - J instead of U - V versus V - J.ln fact, 
NUV - is a better indicator of the current versus past star for- 
mation activity (e.g. Martin et al. 2007, Arnouts et al. 2007fl 
Moreover, the dynamical range covered by the NUV - rest- 
frame colour is larger than the one covered by U-V, making the 
NUV - rest-frame colour less sensitive to unceitainties linked 
to observations. Finally, the NUV rest-frame is still sampled by 
optical data at z > 2 which is no longer true for the rest frame U 
band. 

Figure |3] shows the two-colour criterion used to select the 
quiescent population. Galaxies with MAff/v -Mr > 3{Mr-Mj)+l 
and Mnuv - Mr > 3.1 are considered as quiescent. The advan- 
tage of this classification is that it avoids a mix between dusty 
star-forming galaxies and quiescent galaxies: extinction moves 
star-forming galaxies along a diagonal axes from the bottom left 
to the top right of Figure |3] We derive the absolute magnitudes 
with the method of Ilbert et al. (2005) which minimizes the k- 
correction dependency^. 

The specific star formation rate (sSFR) is the ratio between 
the instantaneous SFR and the stellar mass obtained from the 
best-fit template. In Ilbert et al. (2010) and Dommguez Sanchez 
et al. (2011), we considered that a galaxy was quiescent when 
log(sSFR) < -11. The fraction of galaxies with log{sSFR) < 
- 1 1 in the selection box is given in each panel of Figure |3] Over 
the full redshift range, ~9Q% of the galaxies with log(sS FR) < 
- 1 1 are in the selection box. We also indicate the fraction of 
galaxies with log(sSFR) > -11 in the selection box. The frac- 
tion of galaxies with log{sS FR) > -11 increases from 20% at 
z < 1 to 60% at 2 < z < 2.5. Therefore, both classifications pro- 
vide similar results at z < 1 . But the classification based on the 
sSFR is more conservative at high redshift. 



4. Method used to estimate the Galaxy Stellar Mass 
Function 

We use the tool ALF (Ilbeit et al. 2005) to derive the galaxy stel- 
lar mass functions. This tool includes three non-parametric esti- 
mators (Vmax, SWML, C+) detailed in Appendix A.2 of Ilbert 
et al. (2005). We limit our analysis to above the stellar mass limit 
given in Table 2 (column 4), as defined in ^3.21 We verified that 
the three non-parametric estimators are in good agreement over 
the considered mass ranges. As shown in Ilbert et al. (2004), the 
estimators behave differently if a galaxy population is missing 
and diverge below a given mass limit when a significant pop- 
ulation is no longer present. The good agreement between the 
three estimators confirm that our stellar mass limits are well es- 
tablished. 



4.1. Associated uncertainties 

Our error budget includes the Poissonian errors, the photo- z red- 
shift uncertainties, and the uncertainties in the mass estimations. 

We estimate the cosmic variance with the public tool getcv 
provided by Moster et al. (201 1). The square root of the cosmic 
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Fig. 3. Two-colour selection of the quiescent population. The 
galaxies falling in the black box in the top left are selected as 
quiescent. The fraction (in %) of log(sS FR) < -11 galaxies se- 
lected by the two-colour selection box is the top number in each 
panel. The fraction of galaxies with log{sSFR) > -1 1 within the 
two-colour box is the bottom number. 
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Fig. 4. Fractional error in density as a function of the stellar mass 
and redshift. The top panel shows the errors due to the cosmic 
variance. The middle and bottom panels are the errors associ- 
ated to the template fitting procedure (photo-z and stellar mass) 
for the full sample and the quiescent population, respectively. 
Results are shown only in the mass range covered by our dataset. 



variance for galaxies is computed as the product of the galaxy 
bias and the dark matter root cosmic variance. The dark matter 
cosmic variance is computed for a field of 1 .5 deg^. Moster et 
al. (2011) provide an estimate of the bias as a function of red- 
shift and stellar mass. The fractional error in density is shown 
in the top panel of Figure |4] This uncertainty increases with the 
mass since the massive galaxies are more biased than lower mass 
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z-bin 
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Star-forming 0.20 0.50 22818 7.86 

galaxies 0.50 0.80 26619 8.64 

0.80 1.10 24248 9.04 

1.10 1.50 30058 9.29 

1.50 2.00 19477 9.65 

2.00 2.50 6698 10.00 

2.50 3.00 3060 10.20 

3.00 4.00 1682 10.29 
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Table 2. Best-fit parameters for the full sample, 
forming galaxies. We adopt a double Schechter 
to fit the Vmax data points. A simple Schechter 
star-forming galaxies at z > 1.5. 



for the quiescent galaxies classified with the two-colour criteria and for the star- 
function (see Eq.2) convolved with stellar mass uncertainties (see Appendix A) 
function is considered for the quiescent sample at z > 0.5 and for the global and 



galaxies. The error associated to the cosmic variance is always 
larger than 5% at low mass and reaches 30% at high mass. 

In order to take into account the uncertainties generated by 
the template-fitting procedure (photo-z and stellar mass), we cre- 
ated a set of 30 mock catalogues by perturbing each flux point 
according to its formal error bar. We scatter the photometric red- 
shifts according to their Icr error provided by Le_Phare. For 
each realisation, we recompute the stellar masses, and perform 
again the galaxy classification and recompute the MFs. We com- 
pute the Icr dispersion in density at a given mass and redshift. 
The evolution of this dispersion is shown in Figure |4] These un- 
certainties are below 2% at low masses and z < 2. In this regime, 
the cosmic variance dominates the error budget. At high masses 
Ai > 1O"A10, the errors associated with the template fitting 
procedure reach 30-40% and become as important as the cosmic 
variance. 

In order to obtain the total en^ors, we add in quadrature 
the errors due to the galaxy cosmic variance {errcosmic_var), the 
ones linked to the template-fitting procedure (errfu) and the 
Poissonian errors {errpoissonY 



err,o, 



(1) 



Since errcosmic_var and errf,, are derived as a function of the mass 
and the redshift, it is straightforward to associate the total errors 
to the non-parametric MF data points. 

4.2. Fit of the stellar mass function and Eddington bias 

We fit a parametric form over the Vmax non-parametric data. 
The choice of the non-parametric estimator has no impact on 
our results since we work in a mass range where the three non- 



parametric estimators are consistent. Following Pozzetti et al. 
(2010), we use a double Schechter function, defined as: 



(piM)dM = e af 



M_ 

M* 



M_ 

M* 



dM 



(2) 



with Ai* the characteristic stellar mass, ai and 02 are the 
slopes which need to satisfy the following criteria a2 < ffi, 
and (p* and correspond to the normalization. For the quies- 
cent population, we use a simple Schechter function at z > 0.5 
since we do not detect any upturn at the low-mass end given our 
survey limit. At 0.2 < z < 0.5, we still need to use a fit with a 
double Schechter function. For the global and star-forming pop- 
ulations, we also use a simple Schechter function when the five 
free parameters are no longer constrained (at z > 1.5). At these 
redshifts, a simple Schechter function provides a good fit over 
our considered mass range. 

A crucial step our fitting procedure is to account for the un- 
certainties in the stellar mass. These uncertainties could bias our 
estimate of the high-mass end (Caputi et al. 2011). Since the 
galaxy density exponentially decreases towards massive galax- 
ies, errors in the stellar mass scatters more galaxies into the mas- 
sive end than the reverse (Eddington 1913). Our procedure to 
avoid this bias is detailed in Appendix A. First, we find that the 
stellar mass uncertainties are well characterised by the product 
of a Lorentzian distribution L{x) = ^ (l^+^i with t - 0.04(1 +z) 

and a Gaussian distribution G with cr = 0.5. Then, we convolve 
the double Schechter function (p by the stellar mass uncertain- 
ties: (pconmived - (/> * {L X G). Finally, we fit (pconmived to the Vmax 
non-parametric data. Therefore, the best-fit parameters that we 
provide in Table 2 are deconvolved by the expected stellar mass 
uncertainties and do not suffer from Eddington bias. 
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Fig. 5. Galaxy stellar mass function up to z = 4 for the full sam- 
ple. Each colour corresponds to different redshift bins of vari- 
able step size. Fits are shown in the mass range covered by our 
dataset. The filled areas correspond to the 68% confidence level 
regions, after accounting for Poissonian errors, the cosmic vari- 
ance and the uncertainties created during the template fitting pro- 
cedure. The open triangles and squares correspond to the local 
estimates by Moustakas et al. (2013) and Baldry et al. (2012), 
respectively. 



5. Results: Evolution of the Galaxy Stellar Mass 
Function and Stellar Mass Density 

The galaxy stellar mass functions are computed with a sample of 
220,000 galaxies selected at < 24. We keep only the sources 
in areas with good image quality, representing an area of 1.52 
deg^. We remove the stars and X-ray detected AGNs (Brusa et 
al. 2007). Figure |5] and Figure |6] show the galaxy stellar mass 
functions for the full sample, the quiescent and the star-forming 
populations. The best fit parameters are given in Table 2. In this 
Section, we describe our results out to z = 4. 

5.1. Evolution of the full sample 

Figure |5] shows the evolution of the MFs for the full sample. 
A first option is to consider a pure evolution in stellar mass. 
In this case, we assume that only star formation drives the MF 
evolution (no galaxy can be created or destroyed). We find that 
the evolution is strongly mass-dependent, with the low-mass end 
evolving more rapidly than the high-mass end. For instance, the 
stellar mass of a IO'^^AIq galaxy increases by 0.9 dex between 
1.5 < z < 2 and 0.2 < z < 0.5, while the stellar mass of a 
10^^-^ Mq galaxy increases by only 0.2 dex in the same time 
interval. Therefore, we conclude that the evolution is strongly 
mass-dependent, in agreement Marchesini et al. (2009). A sec- 
ond option is to consider a pure density evolution. A constant 
increase in density by 0.3-0.4 dex, independent of the mass, is 
sufficient to match the 1.5 < z < 2 and the 0.2 < z < 0.5 MFs. 
However such a pure density evolution scenario is not applicable 



star-forming 
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Fig. 6. Galaxy stellar mass function up to z =4 for the star- 
forming population (top panel) and for the quiescent population 
(middle panel). Symbols are the same as Figure |5] The bottom 
panel shows the percentage of quiescent galaxies as a function 
of stellar mass in the same redshift bins. 
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Fig. 7. Stellar mass density as a function of cosmic time. Black 
and red points correspond to the full and quiescent populations, 
respectively. The circles correspond to our new results using 
Ultra VISTA. Solid and open red circles correspond to the two- 
colour and sSFR selected quiescent galaxies, respectively. The 
green shaded area corresponds to the cosmic SFR compiled by 
Behroozi et al. (2012) and integrated over cosmic time assuming 
a return fraction of 0.3. 
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Fig. 8. Stellar mass density as a function of cosmic time. Black 
open circles, red solid circles and blue pentagons correspond to 
the full, quiescent and star-forming galaxies, respectively. 



to the full sample: it would mean that new galaxies which were 
not present in a given redshift appear in the next redshift bin. 
Major mergers are not an option for a pure increase in density 
with cosmic time: for a a = -1.4 MF slope, the density of low 
mass galaxies would decrease by 0.16 dex if we assume that all 
galaxies encounter a major merger since z = 2Q. 

We derive the stellar mass density by integrating the best-fit 
double Schechter functions over the mass range 10** to IO'^AIq. 
The results are shown in Figure Q (black circles). We find that 
the global stellar mass density increases by 1.2 dex between 3 < 
z < 4 and 0.8 < z < 1.1 (a factor 13 in 4.2 Gyr). The evolution is 
much smaller at z < 1 with an increase of 0.2 dex between 0.8 < 
z < 1.1 and the local estimate (a factor 1.6 in 6 Gyr). Therefore, 
the stellar mass is assembled twice as fast at 1 < z < 4 (14% of 
the local stellar mass density per Gyr) than at z < 1 (6% of the 
local stellar mass density per Gyr). Figure |7] shows an excellent 
agreement between the various estimates from the literature, in 
particular at high redshift with the previous results from Caputi 
et al. (201 1), Perez-Gonzalez et al. (2008) and Marchesini et al. 
(2009). We derive stellar mass densities which are lower by 0.2- 
0.3 dex at z > 2 in comparison with Santini et al. (2012). They 
find a slope much steeper than ours {a < -1.8) which explains 
this difference. Their steep slopes are explained by a steep upturn 
in their faintest bins in stellar mass. They use a HAWK-I Kj 
selected sample which is 1 .5 magnitude deeper than our sample 
and we are not able to reach the stellar mass limit at which this 
upturn occurs. Kajisawa et al. (2009) detected a steepening of 
the MF slope with redshift with a reaching -1.6 at 2.5 < z < 3.5. 
While our sample has a similar depth, we find a flatter slope 
with -1.4 < a < -0.9 atz > 1.5, and no clear trend for the slope 
evolution. 

Finally, we compare the stellar mass density with the star for- 
mation rate integrated along cosmic time (green area in Figure 



' The MF would be shifted in density by -0.3 dex (half as many galax- 
ies) and the masses would increase by 0.3 dex 



|7]l. Behroozi et al. (2012) compiled recent estimates of the cos- 
mic SFR. They fit these data with a function (Table 5 of Behroozi 
et al.) that we integrate along cosmic time. We assume that 30% 
of the created stellar mass returns into the interstellar medium 
(Prantzos 2008). We confirm the tension between both measure- 
ments: the integrated SFR is larger by 0.2-0.3 dex than the stellar 
mass density over the full redshift range. 

5.2. Evolution of the quiescent population 

Middle panel of Figure |6] shows the MF evolution of quies- 
cent galaxies. The MF evolution of the quiescent population is 
clearly mass dependent at z < 1 . In this redshift range, we do 
not find any significant evolution of the high-mass end while 
we observe a clear flattening of the faint-end slope. Between 
0.8 < z < 1.1 and 0.2 < z < 0.5, the density of galaxies 
more massive than lO^'^AI© does not increase, while galaxies 
are continuously "quenched" at the low-mass end. For instance, 
the density of IO'^-^Mq galaxies increases by a factor of > 5 
between 0.8 < z < 1.1 and 0.2 < z < 0.5. 

At 1 < z < 3, we find a rapid increase in density of all qui- 
escent galaxies. In contrast with the result at z < 1, the evolution 
is no longer mass dependent. The most massive galaxies evolve 
as fast as intermediate mass galaxies. The density of 10" Al© 
galaxies increases by 1.45 dex (factor 30) between 2.5 < z < 3 
and 0.8 < z < 1.1. The normalisation parameter <!>* increases 
continuously between z ~ 3 and z ~ 1. We consider a pure 
density evolution for the quiescent MF rather than a pure mass 
evolution. Indeed, new quiescent galaxies are created along cos- 
mic time by quenching of star-forming galaxies. Moreover, an 
isolated quiescent galaxy grows by less than 0.03 dex in 6 Gyr 
(since log(sS FR) < -1 1) due to its own star-formation. 

Figure |7] shows the evolution of the stellar mass density for 
the quiescent population. The stellar mass density increases by 
1.5 dex between 2.5 < z < 3 and 0.8 < z < 1.1. Again, we find a 
change of regime around z ~ 1 with the mass assembly slowing 
down at z < 1 . We find that the stellar mass assembly is faster 
at 1 < z < 3 for the quiescent population than for the global 
population which evolves by 0.8 dex in the same redshift range. 
For the massive quiescent galaxies (Ai > lO^AI©), Brammer 
et al. (2011) find an evolution of 0.5 dex between z - 2.1 and 
z = 1. But their mass density drops quickly at z > 2.1. Their 
evolution reaches 1 dex if we consider the redshift range z - 2.3 
and z = 1, which is really close to our value of 1.1 dex. We 
find a slightly slower evolution compared to Ilbert et al. (2010). 
However, our conclusion depends on the method used to classify 
the quiescent population. If we use a classification based on the 
sSFR, as in Ilbert et al. (2010), the stellar mass densities are 
consistent between both studies. The sSFR classification is more 
restrictive than the two-colour selection (see g3.3b . 

5.3. Evolution of the star-forming population 

The top panel of Figure |6] shows the stellar mass function of 
the star-forming galaxies. We can divide this evolution in two 
regimes: above and below 10 '"■^"'" '^Al© (the characteristic stel- 
lar mass, see Table 2). If we consider only the low mass regime, 
the faint-end slope remains steep over the full redshift range. The 
slopes stay in the range -1.5 < c < -1.3 except at 2 < z < 2.5 
where the slope is shallower We do not detect any significant 
trend in the slope evolution over the stellar mass range covered 
by our dataset. We observe a strong evolution of low-mass galax- 
ies, especially at z > 2. For the MF at 2 < z < 2.5 to match that at 
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Fig. 9. Evolution of the sSFR (SFR/M). The black filled cir- 
cles correspond to the sSFR derived from the evolution of the 
star-forming MF evolution for low-mass galaxies, and should be 
considered as a lower limit. The other points are a compilation 
from previous studies. For Karim et al. (201 1), we used the sSFR 
measured in the mass bin 10'"'^ - 10"'-^ Al©. 

1.5 < z < 2 requires the stellar mass of a lO'^'^At© galaxy to in- 
crease by around 0.4 dex. We discuss in ^6.1l how such evolution 
could be interpreted in term of sSFR. If we consider now the evo- 
lution of the high-mass end, we do not detect an evolution of the 
density of the most massive galaxies 10^ ^■^^^^'^Mq. Since these 
star-forming galaxies are forming new stellar populations, these 
massive star-forming galaxies are necessarily quenched along 
cosmic time, as we will discuss in 36.21 

Figure [8] shows the evolution of the stellar mass density for 
star-forming galaxies. As for the quiescent population, we also 
observe two regimes with a faster evolution at 1 < z < 4 than at 
z < 1. We observe an evolution of 0.5 dex between 2.5 < z < 3 
and 0.8 < z < 1 . 1 , while the quiescent galaxies evolve by 1 .5 dex 
in the same period. Therefore, quiescent galaxies are building 
faster at this epoch. We also observe that the stellar mass density 
of quiescent and star-forming galaxies are comparable at z < 1, 
while star-forming galaxies dominate the stellar mass budget at 
higher redshift. 

6. Discussion 

6.1. Linking the star-forming MF evolution witli the sSFR 
evolution 

We now consider an admittedly over-simplistic scenario in 
which evolution of star-forming galaxies is driven only by star 
formation (i.e. we consider that mergers do not significantly 
change the galaxy distribution between two redshift bins, see 
Section 4.3 of Boissier et al. 2010). Given this assumption, 
the star formation increases the stellar mass by Alog(A^) - 
log(l + sSFR(z, M) X dt) in the time interval 6t. Given the tight 
correlation between the SFR and the stellar mass (e.g. Noeske 
et al. 2007, Daddi et al. 2007, Elbaz et al. 2007), we add the 



assumption that all galaxies have the same sSFR at a given red- 
shift. This assumption is not valid at M > 1O'°'^"'°'^A10 since 
we show that the MF evolution is mass dependent and the sSFR 
is reduced for massive galaxies (Karim et al. 2011). However, 
we can still make this assumption for low-mass galaxies. This 
is possible since the slopes are consistent. Given these assump- 
tions, the MF is shifted equally along the mass axes by Alog(AI). 
Let us consider the MF(fi) measured at zi < z < Z2 and the 
MF(f2) measured at Z2 < z < Z3, with zi < Z2 < Z3. We compute 
the shift Alog( Al) which needs to be applied in order to superim- 
pose MF(f2) over MF(fi) at the reference mass M = 10'°-^ Al©. 
With this simple exercise, we convert Alog(A4) into sSFR in 
each redshift bin. This sSFR estimate is a lower limit. Indeed, 
star-forming galaxies could be quenched and move to the quies- 
cent population in the time interval 6t. In this case, Alog(AI) is 
underestimated. The sSFR computed with this indirect method 
(black filled circles in Figure |9]l is compared with direct mea- 
surement of the sSFR from the literature. At z > 2, the contribu- 
tion of the quiescent to the global population is lower than 20%. 
Therefore, we do not expect to underestimate significantly the 
sSFR. At z > 2, we find that the sSFR is in good agreement with 
Karim et al. (201 1), Reddy et al. (2012) and Stark et al. (2012). 
We conclude that the evolution of the low-mass galaxies that 
we observe could be fully explained by a star-formation activ- 
ity consistent with the sSFR derived in the literature. At z < 2, 
the quiescent galaxies are built very fast and their contribution 
to the global MF reaches 50%, which explains why the sSFR we 
get are so lower in comparison to the literature. 

6.2. Quenching processes 

Newly quenched galaxies and major mergers drive the evolution 
in density of the quiescent MF, since by definition a quiescent 
galaxy grows by less than 0.03 dex in 6 Gyr by star formation. 
Therefore, the efficiency of quenching with cosmic time can be 
directly estimated from the MF evolution of the quiescent and 
star-forming galaxies. In this section we discuss our results in 
term of quenching efficiency, and we summarize our interpreta- 
tion with a synthetic scheme (FigurefTOli. 

As we discuss in 36.11 the evolution of the low-mass star- 
forming galaxies at z > 2 is consistent with the expected sSFR 
while the lack of evolution of the massive end implies that mas- 
sive star-forming galaxies are necessarily quenched. Our data are 
deep enough at z < 3 to show that the slope of the star-forming 
sample is steeper than the slope of the quiescent population (see 
Figure |6]l. If the quenching processes acted with the same effi- 
ciency at all masses, the MF of the quiescent population should 
have the same slope as the star-forming MF. We conclude that 
the physical process which quenches star formation is necessar- 
ily more efficient above M > lO'^'^^'^'^A^©, i.e. the maximum 
in density of the quiescent MlQ. This scenario is fully consistent 
with the model proposed in Peng et al. (2011) who introduce a 
"mass quenching" process which is independent of environment 
for the massive galaxies. We show that this mass quenching has 
been active since z = 3. This mass quenching process must have 
been efficient at 1 < z < 3 to rapidly build up the quiescent popu- 
lation (left panel of Figure [TOli. For instance, around 30% of the 
star-forming galaxies with Ai ~ lO^'AI© should be quenched 
between 2.5 < z < 3 and 2 < z < 2.5 to explain a density 
increase by 0.25 dex of massive quiescent galaxies. While we 
can see that this quenching process is extremely efficient above a 
given mass, we cannot infer the physical nature of the quenching 
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Fig. 10. A possible scenario showing how sSFR and quenching impact the star-forming MF (blue lines) and quiescent MP (red 
lines). The large arrows represent quenching. The cyan arrows correspond to the mass increase expected in 2 Gyr, by taking the 
sSFR values compiled by Stark et al. (2012). The left panel corresponds to the high redshift bin 2.5 < z < 3, where we show the 
mass quenching process which is more efficient at high mass. The middle panel corresponds to the redshift bin 1.1 < z < 1.5 when 
the quiescent starts to dominate the high-mass end and the reservoir of massive star-forming galaxies become depleted. The right 
panel corresponds to the redshift bin 0.2 < z < 0.5, when environment quenching generates new low mass quiescent galaxies. 



process. AGN feedback or major mergers (or both) are possible 
candidates (e.g. Hopkins et al. 2006). 

Since the star-forming galaxies dominate the high-mass end 
of the global MF at z > 1 (bottom panel of Figure|6]l, an efficient 
"mass quenching" process is rapidly transitioning them into pas- 
sive systems. If the "mass quenching" process depletes more 
rapidly the reservoir of massive star-forming galaxies than new 
star formation is able to replenish it, we reach an epoch where 
the quenching is no longer a channel to create new massive qui- 
escent galaxies. By z ~ 1 the supply of massive star forming 
galaxies has dwindled to the point that few new galaxies are be- 
ing quenched. It does not mean that the quenching process is not 
efficient anymore, but there is no massive star forming galaxies 
to be quenched. Figure |5] shows that > 70% of the galaxies more 
massive than M > 10" Al© are quiescent at z < 1.1. Moreover, 
the sSFR estimated from the literature decreases by a factor 15 
between z ~ 1.5 andz ~ 0.6 (see Stark et al. 2012, Figure|9)- We 
conclude that the reservoir of massive star-forming galaxies has 
been depleted at z < 1 and that the growth rate of star-forming 
galaxies is no longer sufficiently efficient to generate numerous 
massive star-forming galaxies (middle panel of Figure [TOt. 

The flattening of the faint-end slope of the quiescent MF at 
z < 1 (see ^5.2l i shows that the quenching starts to be efficient at 
low masses (right panel of FigurefTOt. Using the COSMOS sam- 
ple, Scoville et al. (2013) find that low-mass quiescent galax- 
ies appear in high density environment at z < 1. Peng et al. 
(2011) proposed the "environment quenching" process which 
mainly affects the satellite galaxies as large-scale structure de- 
velop. Therefore, we interpret the flattening of the slope of the 



quiescent MF as an environmentally driven process, even if our 
results alone can not demonstrate it. 



6.3. Impact of the assumed star formation histories 

We derive the stellar masses assuming an exponentially declin- 
ing star formation rate for the stellar population synthesis mod- 
els. However, numerous studies pointed out that such star forma- 
tion histories could be inappropriate at z > 1 (e.g. Maraston et 
al. 2010, Pacifici et al. 2013). Therefore, we investigate how this 
choice impacts our conclusions. 

We derive the MF using two other libraries: one library based 
on delayed star formation histories and another using a compos- 
ite of various star formation histories (see ^3.2l l. We generate the 
second library based on delayed star formation histories follow- 
ing a law in r^^fe"'^'". The maximum Star Formation Rate (here- 
after SFR) peak could be delayed by 0.1, 0.5, 1, 2 or 3 Gyr and 
we consider three different metallicities (Z = 0.004, Z - 0.008, 
Z = 0.02). We generate a third composite library including vari- 
ous combinations: three templates with an exponentially declin- 
ing SFR (t = 0.1, 1,5 Gyr), two templates with a delayed star 
formation history (maximum peak at 1 and 3 Gyr). For each of 
these templates, we allow a possible second burst 3 Gyr after 
the peak in SFR. The amplitude of this second burst could be 
0.25 or 0.5 times the maximum of the first burst. This second 
burst is generated using a constant SFR. We derive the MFs with 
these three libraries (called hereafter exponential, delayed and 
composite libraries). 
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Figures[TT]and[T2]show the MFs of the star-forming and qui- 
escent galaxies, respectively. The MFs based on the exponential, 
delayed and composite libraries are shown with red squares, blue 
triangles and green crosses, respectively. We do not detect any 
significant difference between the MFs. Therefore, all our con- 
clusions are equally valid for any of these three libraries. 

In detail, we find that quiescent galaxies are almost always 
fitted with a delayed SFR history peaking after 0. 1 or 0.5 Gyr, 
which is close to an exponentially declining history. For the star 
forming galaxies, the templates are uniformly distributed be- 
tween the five delayed star formation histories. By comparing 
the best-fit ages and the t values of the delayed templates, 30% 
of star-forming galaxies at z > 2 are fitted with a template in the 
rising part of the star formation history. The dispersion between 
the stellar masses computed with the exponential and delayed li- 
braries is lower than 0.06 dex out to z < 4, with no systematic 
difference. This explains why the delayed and exponential MFs 
are similar For the composite library, the templates with an ex- 
ponentially declining star formation history are selected in 87% 
of the case, with a possible small second burst. This explains 
why we obtain almost the same results between the composite 
and the exponential libraries. 

6.4. Comparison with the semi-analytical models 

We now compare our mass function measurements with the pre- 
dictions of semi-analytical models. The mock catalogs are based 
on ACDM simulations from Wang et al. (2008) with the cos- 
mological parameters derived from the third-year WMAP data 
{Ho = 74.3 km/s, Qm = 0.226 and Qa = 0.774). The Hght cone 
survey covers an area of 1.4 X 1.4 deg2 similar to COSMOS. 
Galaxy properties were generated using the galaxy formation 
model, as detailed in De Lucia & Blaizot (2007) and Wang & 
White (2008). Since the redshift and the galaxy stellar masses 
are available for all galaxies in the simulation, we can directly 
compute the predicted MFs. We use our standard cosmology 
(Hq = 70 km/s, Q.M - 0.3 and Sl\ - 0.7) to renormalise the 
predicted stellar mass counts by the comoving volume. We iso- 
late the quiescent population in the simulations. We try three 
different criteria: a red clump criteria using rest-frame UV - r 
versus r - K colours, a cut in log{sS FR) < -11 and a cut in 
B - I corrected for dust extinction. The quiescent MFs are al- 
most identical regardless of selection criteria. 

The comparison between the observed and predicted MFs is 
shown in Figure [TT] for the star-forming sample. First, we find 
that the predicted and observed faint-end slopes are in good 
agreement at least up to z - 2. Secondly, the model under- 
predicts the density of massive galaxies. The mismatch between 
the observed and predicted high-mass end increases with red- 
shift. If we add a Icr error of 0.2 dex to the stellar masses, as 
commonly assumed in the literature (e.g. Cattaneo et al. 2006, 
Bower et al. 2012), we reduce the tension between model and 
observation (dotted lines in Figure fTTT i. We point out that the 
0.2 dex uncertainty commonly associated to the stellar mass 
is mostly explained by systematic uncertainties (choice of the 
IMF, adopted population synthesis models) and should not be 
applied to the stellar masses. For a given library, the statistical 
errors are much smaller Moreover, we already deconvolved the 
fit by the expected statistical errors in the stellar masses (see 
®. Therefore, the red solid line is directly comparable to the 
model prediction, showing that the predicted high-mass end is 
underestimated. Several factors could explain this mismatch: 1) 
our high redshift sample is contaminated by low redshift galax- 
ies (i.e. catastrophic failures in the photo-z); however Figure [1] 



do not show such contamination with our current spectroscopic 
sample; 2) some assumed physical laws (e.g. IMF) evolve with 
redshift which generate redshift dependent biases in our stellar 
mass estimates; 3) the models do not generate enough massive 
galaxies at z = 3 - 4. The last option would require that the mod- 
els generate 10"^ Al© galaxies in less than two Gyr Then, these 
galaxies should stop to grow on a very short timescale by having 
their star formation activity drastically reduced or quenched. 

Figure [12] shows the comparison between observed and pre- 
dicted MFs for the quiescent population. The mismatch be- 
tween the observed and predicted faint end slopes is dramatic. 
At 9.5 < log(M) < 10 and 0.5 < z < 0.8, the model over- 
estimates the low-mass end by a factor 10. Such an effect has 
been already noticed by Cucciati et al. (2012) in the VVDS field 
based on the B-band luminosity function. Wang & White (2008), 
Bielby et al. (2012), Bower et al. (2012), Guo et al. (201 1) show 
that the models overproduce the density of low-mass galaxies 
for the full population. In this paper, we show that the origin 
of this discrepancy is the high density of low-mass quiescent 
galaxies, since the faint end slope of the star-forming galaxies is 
correctly predicted at z < 2. Information such as the host halo 
mass, the galaxy location in the halo (central or satellite) and 
the galaxy morphology are kept in the simulation. We use this 
information to investigate in more detail the properties of the 
quiescent galaxies overproduced by the models. We first keep 
only bulge-dominated quiescent galaxies with a bulge-to-disk 
ratio greater than 0.5 (brown short dashed lines). The slopes are 
in much better agreement with the data. Therefore, the model is 
able to generate a reasonable amount of quiescent galaxies when 
the quenching is associated to a process which makes the galaxy 
dominated by a bulge. Indeed, Ilbert et al. (2010) showed that 
between 40 and 70% of the quiescent galaxies have an ellipti- 
cal morphology at z < 0.8 and 9.5 < log(Ai) < 10, while in 
the model, 20% of the quiescent galaxies are bulge dominated 
in the same mass range. We also show the predicted MF of cen- 
tral quiescent galaxies (cyan long dashed lines). It appears that 
the low-mass quiescent simulated galaxies are mostly dominated 
by satellite/orphans galaxies. One possible interpretation of the 
mismatch between the faint-end slopes is an over-quenching of 
the star formation in satellite galaxies and/or an over-quenching 
of the star-formation in disk galaxies. 



7. Conclusions 

We study the stellar mass assembly out to z = 4 with a unique 
sample of 220,000 galaxies selected at K, < 24 in the COSMOS 
field. Galaxies are selected using the new UltraVISTA DRl near- 
infrared data release, which reaches between one and two magni- 
tudes deeper than previous COSMOS near-infrared data. Fluxes 
are measured on PSF-homogenised images in 25 bands includ- 
ing optical broad and intermediate-band Subaru data taken for 
the COSMOS project. We also add four bands at 3.6 - Spm 
from S-COSMOS and the GALEX NUV band. We keep only 
the sources at < 24 located in areas with a good image qual- 
ity, which represents an effective area of 1.52 deg^. 

The photometric redshifts are derived using Le_Phare fol- 
lowing a procedure similar to Ilbert et al. (2009). Our photo- 
metric redshifts have two regimes: at i^^ < 22.5 {Zmed ~ 0.5), 
precision 1% with less than 1% of catastrophic failures. In the 
high redshift range 1.5 < z < 3.5, the precision of the photo- 
z is tested against the zCOSMOS faint sample, and is 3% for 
Cicrf ~ galaxies. We estimate the galaxy stellar masses using 
a library of synthetic spectra generated using the SPS model of 
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Fig. 11. Galaxy stellar mass functions of the star-forming sample. The points correspond to the Vmax estimator (we do not show the 
other estimators for clarity) and error bars include cosmic variance, Poissonian errors and photo-z/stellar mass uncertainties. The 
red squares, the blue triangles and the green crosses correspond to the MFs computed using the exponential, delayed and composite 
libraries. The red solid line is the best-fitted double-Schechter function. The black solid, dashed and dotted lines correspond to the 
MFs predicted by the semi-analytical model, assuming no error on the stellar masses, an error evolving as described in appendix A 
(using the product of a Lorentzian and gaussian distribution), and a gaussian with a cr = Q.ldex. 



Bruzual and Chariot (2003). We show that our results are inde- 
pendent of the assumed star-formation history. 

We find that the evolution of the global stellar mass func- 
tion is strongly mass-dependent. The low-mass end evolves more 
rapidly than the high-mass end. If we consider an evolution 
purely driven by star formation, the low-mass galaxies evolve 
by almost 1 dex between 1.5 < z < 2 and 0.2 < z < 0.5, while 
the stellar mass of the most massive galaxies increases by less 
than 0. 1 dex in the same time interval. The observed evolution 
for low-mass galaxies is fully consistent with the expected sSFR 
measured at z > 2. The lack of evolution of the massive end can 
be interpreted as a direct consequence of star formation being 
drastically reduced or quenched when a galaxy becomes more 
massive than M > 10'°-^""' '^Mo. 



We derive the comoving stellar mass density of the full sam- 
ple from z = 4 to z = 0.2. We find that the stellar mass is as- 
sembled twice as fast at 1 < z < 4 (14% of the local stellar 
mass density per Gyr) than at z < 1 (6% of the local stellar mass 
density per Gyr). We also confirm that the SFR integrated along 
cosmic time overproduces by 0.2-0.3 dex the stellar mass density 
at all redshifts, even using the latest cosmic SFR compilation by 
Behroozi et al. (2012). 

We isolate the quiescent population using a classification 
based on the rest-frame colours NUV - and - J. This clas- 
sification separates cleanly dusty star-forming galaxies and qui- 
escent galaxies. We find that the MF evolution of the quiescent 
population is definitively mass-dependent at z < 1, confirming 
the trend shown in Ilbert et al. (2010). We do not find any sig- 
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Fig. 12. Galaxy stellar mass functions of the quiescent sample for selected using a two-colour technique. Symbols are the same 
as Figure [TT] The short dashed brown lines correspond to quiescent galaxies with a bulge-to-disk ratio greater than 0.5. The long 
dashed cyan lines correspond to the central galaxies in the models. 



nificant evolution of the high-mass end at z < 1 (any evolution 
being limited to AM < Q.ldex), while we observe a clear flat- 
tening of the slope in the same redshift range. We interpret this 
evolution of the low-mass end of the MF as arising from contin- 
uous quenching of galaxies between z ~ I and z ~ 0.1, probably 
by physical processes related to environment. 

The Ultra VISTA data now allows us to trace the growth in 
stellar mass density in quiescent galaxies from z ~ 3 to the 
present day. From z ~ 3 to z ~ 1 we find a rapid increase in 
the stellar mass density of all quiescent galaxies, independent 
of stellar mass. We confirm that the steep rise of more than 
one order of magnitude of this population between 1 < z < 2 
observed in previous works (Amouts et al., 2007; Ilbert et al., 
2010) follows an earlier rise by a factor ~3in2<z<3.5. 
In total, we find that the density of 1O"A10 galaxies is contin- 
uously increasing, reaching a factor ~ 30 between 2.5 < z < 3 
and 0.8 < z < 1.1. This must indicate that a fraction of star- 



forming galaxies is continuously quenched at z > 1. Because 
of the different faint-end slopes of the MF between the quies- 
cent and global populations at z < 3, we infer that the physi- 
cal processes which quench the star formation are more efficient 
above Ai > 10'"'^"'°'^A1o (the maximum in density of the qui- 
escent MF). This scenario is consistent with the model proposed 
by Peng et al. (2011) who introduce a "mass quenching" pro- 
cess. Since the high-mass end of the quiescent MF stops evolv- 
ing at z < 1, we conclude that: 1) star formation is not efficient 
enough at z < 1 to produce new massive star-forming galaxies, 
which could be quenched later; 2) major mergers between mas- 
sive galaxies are not sufficiently frequent at z < 1 to increase 
significantly the density of massive quiescent galaxies. 

Finally, we compare our data with the predictions of the 
semi-analytical model of De Lucia and Blaizot (2007) as de- 
scribed in Wang et al. (2008). We find that the model under- 
predicts the density of massive galaxies for quiescent and star- 
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forming galaxies. The tension between model and observations 
at the high-mass end can be alleviated if scatter in the modeled 
stellar masses errors is a Icr of 0.2 dex, which is larger than the 
expected noise in the stellar mass estimate. When we compare 
the observed and predicted MFs for the quiescent population, the 
difference in the slope is striking. At logiM) ~ 10 and z < 0.8, 
the model overestimates the low-mass end by a factor 10. On 
the other hand, the model successfully reproduce the faint-end 
slope of the star-forming population. This indicates that recipes 
to quench satellite and/or disk galaxies within semi-analytical 
models are stUl incomplete and will need to be modified. 

The COSMOS field will be observed during five years for 
the Ultra VISTA survey, allowing to obtain NIR data two magni- 
tudes deeper. In the meantime, the SPLASH program (Capak et 
al., in preparation) will observe the COSMOS field with IRAC, 
allowing a gain of 1.5 mag at 3.6//W. The combination of this 
new data will allowed to extend this work at 4 < z < 7 and better 
constrain the contribution of low-mass galaxies. 
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Appendix A: Treatment of the Eddington bias 

Caputi et al. (2011) outlined that a bias (Eddington 1913) could 
aff'ect the estimate of the high-mass end of the stellar mass func- 
tion. The uncertainties on the stellar masses scatter the galaxies 
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from one stellar mass bin to another In the case of the stellar 
mass function, the density exponentially decreases toward mas- 
sive galaxies. The stellar mass scattering therefore moves more 
galaxies into the massive end than the reverse. As a consequence, 
we overestimate the density of massive galaxies. 

Caputi et al. (201 1) estimated that the Eddington bias could 
affect the density by up to 0.13 dex for the most massive galax- 
ies. In this paper, we go one step further and we deconvolve the 
density estimate by the stellar mass uncertainties. The first step is 
to characterise the stellar mass uncertainties. We do not consider 
here systematic errors linked to the choice of the stellar syn- 
thesis models or to the IMF. We consider only the uncertainties 
linked to the photometry and to the photometric redshifts. We 
generate a "noisy" sample by adding noise to the apparent mag- 
nitudes and to the photometric redshifts using their associated er- 
ror bars. Then, we recompute the stellar masses with this "noisy" 
catalogue. Figure lATI shows the difference between the original 
masses and the ones obtained on this "noisy"catalogue. The red 
histogram shows the cumulative distribution of these differences. 
We first model the stellar mass uncertainties using a Gaussian 
distribution G(x) - — i=ex/7"^7 as it is commonly assumed 
(e.g. Cattaneo et al. 2006, Bower et al. 2012). The dotted lines 
correspond to a Gaussian distribution with cr - 0.04(1 + z). The 
fraction of galaxies with a difference larger than 0. 1 dex is un- 
derestimated and the Gaussian distribution is not a good repre- 
sentation of the main peak at z > 2.5. We also model the distri- 
bution by using a Lorentzian distribution L(x) = ^ (l^+k- ^^'■^ 
T = 0.04(1 H-z) (short dashed lines in Figure lATTT l. This parametri- 
sation of the stellar mass uncertainties produces too many galax- 
ies with differences larger than 0.2 dex. Finally, we adopt a func- 
tion which is the product between a Lorentzian distribution with 
T = 0.04(1 -Hz) and a gaussian distribution with cr = 0.5. This dis- 
tribution is shown with the long dashed lines in Figure lATTT l and 
provides the best representation of the stellar mass uncertainties. 

The second step is to take into account the stellar mass un- 
certainties when we derive the stellar mass function. We decon- 
volve only the parametric fit of the stellar mass function by the 
stellar mass uncertainties (not the non-parametric estimates). As 
described in ^ we fit a double-Schechter function (piAi) to the 
Vmax non-parametric data. In order to take into account the stel- 
lar mass uncertainties, we convolve the double Schechter func- 
tion by our modelized stellar mass uncertainties 

^ix)L{M - X, t)G(M - X, o-)dx (A.l) 

oo 

with T - 0.04(1 + z) and cr - 0.5. We fit the convolved double 
Schechter function to the Vmax non parametric data. Figure [A2] 
shows the result of this convolution at z > 0.8 for the global 
sample. The thick solid black lines show the best fit double- 
Schechter function, deconvolved by the stellar mass uncertain- 
ties (our fiducial measurement used in this paper). The differ- 
ence between the best-fit double-Schechter functions with and 
without the deconvolution (black and green lines, respectively) 
is significant only at the high-mass end and remains below 0. 1 
dex. The blue solid dashed lines correspond to the fiducial MF 
convolved by the G x L functions. By construction, the blue 
dashed line provides a good fit of the Vmax non parametric data 
points. The pink dotted lines correspond to the fiducial MF con- 
volved by a simple gaussian with cr = 0.04(1 + z). As seen in 
the 0.8 < z < 1.1 redshift bin, a convolution by a gaussian does 
not explain the presence of some massive galaxies, which is well 
explained if we use the GxL convolution functions. If we allow 
more extreme errors in the stellar mass estimate, using for in- 



gauss lorentz lorentz*gauss 




log(original masses) - log(noisy masses) 

Fig. A.l. The red histogram show the cumulative distribution of 
the difference between the original stellar masses and the ones 
obtained with a "noisy" catalogues (redshifts and magnitudes are 
scattered according to the expected errors). The dotted, short- 
dashed and long-dashed lines correspond to a Gaussian distri- 
bution with cr = 0.04(1 + z), a Lorentzian distribution with 
T - 0.04(1 + z), and the product of a Lorentzian and a Gaussian 
distributions with r = 0.04(1 + z) and cr - 0.5. 

stance a convolution by the same Lorentzian function multiplied 
by a Gaussian with cr = 1, we would observe a larger population 
of massive galaxies (red long dashed line). 



16 



O. Ilbert et al.: Mass assembly in quiescent and star-forming galaxies since z - 4 from Ultra VISTA 



-3 
-4 

X 



-5 
T3 -6 



-5 



-6 



Ll_lgM ' 1 1 1 1 1 1 1 1 1 1 1 : 

r 

: 0.8-1.1 Y |i! 

' 1 1 1 1 1 1 1 1 1 i\ M ; 


- 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 . 
: 1.1-1.5 

: 1 1 1 1 1 1 1 1 1 ; 


: 1.5-2.0 

: 1 1 1 1 1 1 1 1 1 1 ' 1 1 1 ': 


: 

: 2.0-2.5 v.. \ — : 
— \ ■. ^ — 

: 1 1 1 1 1 1 1 1 i\' 1 : 


- v - 

_2. 5-3.0 \ 


^3.0-4.0 1 Vj)^-^ 



log(M,) [Mq ] 

Fig.A.2. Global stellar mass function in several redshift bins. 
The open circles correspond to the Vmax non-parametric es- 
timates. The thick black solid lines correspond to the double- 
Schechter fit deconcolved by the stellar mass uncertainties. The 
green thin solid lines correspond to the simple fit without de- 
convolution. The pink dotted lines, the short-dashed blue lines 
and the long dashed red lines correspond to the fiducial mass 
function (black solid lines) convolved with G(cr - 0.04(1 -i- z)), 
L(t = 0.04(1 + z))G(o- = 0.5), L(t = 0.04(1 + z))G(o- = 1), 
respectively. 



